#ifndef GENERATOR_DIRICHLET_HEADER
#define GENERATOR_DIRICHLET_HEADER

#include<string>
#include<fstream>
#include<sstream>
#include<iostream>
#include<algorithm>
#include<iterator>
#include<vector>
#include<map>
#include<cstdlib>
#include<ctime>
#include<cmath>
#include <gsl/gsl_randist.h>
#include "def.hpp"

// gsl_ran_dirichlet (const gsl_rng * r, size_t K, const double alpha[], double theta[])

using namespace std;

int GE_Dirichlet(gsl_rng* r){

  ofstream outFile;
  outFile.open(PuzaDEF::Instance()->output_file_name.c_str());

  // recover dirichlet parameters from base measure and scalar
  double* temp_alpha = new double[PuzaDEF::Instance()->GE_Dirichlet_K];

  for(int k=0; k < PuzaDEF::Instance()->GE_Dirichlet_K; k++){
    temp_alpha[k] = PuzaDEF::Instance()->GE_Dirichlet_Alpha * PuzaDEF::Instance()->GE_Dirichlet_Base[k];
  }

  for(int i =0; i < PuzaDEF::Instance()->GE_Samples; i++){
    gsl_ran_dirichlet(r, PuzaDEF::Instance()->GE_Dirichlet_K, temp_alpha, PuzaDEF::Instance()->GE_Dirichlet_Data);
    outFile << i << "\t";
    for(int k=0; k < PuzaDEF::Instance()->GE_Dirichlet_K; k++){
      outFile << k << ":" << PuzaDEF::Instance()->GE_Dirichlet_Data[k] << " ";
    }
    outFile << endl;
  }
  outFile.close();

  delete[] temp_alpha;
  return 0;
}

#endif
